function XZsup=sup(Xsup, Zsup)
[cXsup, cZsup] = ndgrid(Xsup,Zsup);
XZsup = [cXsup(:), cZsup(:)]; %support of (X,Z)  (size(Xsup,2)*size(Zsup,2)x1
end